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ABSTRACT 

We investigate the kinetics of diffusion-controlled heterogeneous single-species annihilation, 
where the diffusivity of each particle may be different. The concentration of the species 
with the smallest diffusion coefficient has the same time dependence as in homogeneous 
single-species annihilation, A + A^O. However, the concentrations of more mobile species 
decay as power laws in time, but with non-universal exponents that depend on the ratios 
of the corresponding diffusivities to that of the least mobile species. We determine these 
exponents both in a mean-field approximation, which should be valid for spatial dimension 
d > 2, and in a phenomenological Smoluchowski theory which is applicable in d < 2. Our 
theoretical predictions compare well with both Monte Carlo simulations and with time 
series expansions. 

P. A. C. S. Numbers: 02.50.-r, 05.40.+j, 82.20.-w, 82.20.Wt 
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I. Introduction 



The kinetics of diffusion-controlled single-species annihilation, A + A — > 0, when 
each particle has the same diffusion coefficient is now well understood [1]. For spatial 
dimension d > 2, the kinetics may be accounted for by the rate equation which predicts 
that the density decays as in the long-time limit. For d < 2, various phenomenological 
approaches predict that the density decays as t"*^/^, but with logarithmic corrections 
appearing in d = 2. Accompanying this relatively slow kinetics is a spatial organization 
in which the probability of finding particles at small separations is reduced compared to 
a random distribution. In one dimension, exact solutions, either based on an occupation 
number formalism [2], or by mapping the reaction onto the kinetic Ising- Glauber model 
[3], provide definitive results about this spatial organization and the reaction kinetics. 

Our goal in this paper is to describe the kinetics of heterogeneous single-species an- 
nihilation, which is defined by the reaction scheme Ai + Aj. 0. Here Ai denotes the 
i^^ species and the reaction rate matrix Kij is a function of the diffusivities of the two 
reacting species. Although we refer to different species in the context of their diffusivity, 
the reaction itself is single-species annihilation with distinct rates for different reaction 
channels. We will consider reactant diffusivities which are drawn from a probability dis- 
tribution. Such a situation arises naturally when the reactants have different masses. As 
we shall show, the kinetics of the heterogeneous system is considerably richer than that of 
the homogeneous analog. Similar behavior was also encountered in heterogeneous single- 
species annihilation with ballistic particle motion [4], where the kinetics depends in an 
essential way on the form of the initial distribution of velocities. The present investigation 
is a natural counterpart of this earlier work for diffusive single-species annihilation. 

When the number of species is finite {i. e., the density distribution, P{D,t), contains 
a finite number of discrete peaks), the rate equations predicts that the least mobile species 
(with diffusivity -Dmin) decays as t~^, as in homogeneous annihilation. However, the more 
mobile species each decay at a faster power-law rate with an associated exponent that 
depends on the diffusivity ratio between the more mobile and the slowest species. When 
P{D,t) is continuous but with -Dmin > 0, the rate equations again show that the least 
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mobile species predominates in the long time limit, and that the decay of the more mobile 
species is described by non-universal power-law behavior. However, the detailed form 
of P{D, t — 0) near -Dmin contributes to a logarithmic prefactor in the decay law. If 
-Dmin = 0, the kinetics strongly depends on the initial conditions. For an initial distribution 
of diffusivities with a power law tail, P{D,0) ~ D'*, as D — > 0, the concentration and 
the average diffusion coefficient decay as c ~ t~°' and (D) ~ t~^, respectively, with a = 
(2 + 2|u)/(3 + 2/x) and (3 — 1 — a — 1/(3 + 2/x). All of these predictions are expected to 
apply for (i > 2. 

For d < 2, we apply the Smoluchowski theory to account for the reaction kinetics. 
We first test this phenomenological approach on the "impurity" problem, in which there 
is a small concentration of particles with one diffusion coefficient in a homogeneous back- 
ground of particles with a different diffusivity. The Smoluchowski theory predicts that the 
impurity species decays as a non-universal power-law in time. These predictions are in 
good agreement with numerical results that arise from time-series expansions and Monte 
Carlo simulations in one dimension. Connections are also made with the exactly-soluble 
problem of the decay of impurities in the related single-species coalescence process. We 
then generalize the Smoluchowski approach to both continuous distributions of diffusivities 
and general spatial dimension d <2. 

The rest of this paper is organized as follows. In Section 2, we derive basic results from 
the rate equations for both discrete and continuous distributions of diffusivities. In Section 
3, we investigate the reaction kinetics for d < 2 by the Smoluchowski approach. This leads 
to phenomenological rate equations with time-dependent reaction rates. The predictions of 
this approach are then outlined. In Section 4, we present time series expansion and Monte 
Carlo simulation results in one dimension to support our theoretical findings. Finally, we 
give a brief discussion in Section 5. 

II. Rate Equations Approach 

Consider first the mean-field rate equations for the case of two distinct species, Ai 
and ^2, with respective diffusivities Di and D2. The two species interact according to the 
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bimolecular processes A^+Ai — > 0, Ai + A2 ^ 0, and A2 + A2 ^ 0, with respective rates 
Kii, K12 = ^^21 5 and K22- The corresponding rate equations are 

^ = -Kiic\ - K12C1C2, (la) 

^ = -K22CI - K12C1C2, (lb) 

where Cj denotes the concentration of the i^^ species. 

Smoluchowski suggested a simple way to relate the reaction rate Kij to the diffusion 
coefficients Di and Dj, and radii Ri and Rj of the reactants [5]. This derivation is based 
on considering the low-density limit. In the rest frame of a particle of species z, it may 
be considered as a spherical stationary trap of radius R which is surrounded by cloud of 
j particles which are captured upon contact with the trap. The reaction rate is identified 
as the flux of particles of type j to the trap under the boundary conditions of absorption 
at the surface of the i^^ particle and a fixed concentration as r ^ 00. Upon solving the 
diffusion equation in three dimensions under these conditions, one straightforwardly finds 
that this flux equals K = AttDR in the long-time limit (see, e. g., [6]). If both species 
perform independent Brownian motions with respective diffusion coefficients Di and Dj, 
then the appropriate generalization of the reaction rate is K^j = 4:7r{Di + Dj){Ri + Rj). 
Assuming, for simplicity, that Ri = R2 = R and absorbing the numerical factor 8ttR into 
the overall time scale, we may re- write Eq. (1) as, 

= -2Dicl - (Di + D2)ciC2, (2a) 

^ = -2D2cl - (Di + D2)ciC2. (26) 

To solve these equations, consider first the time dependence of the concentration 
ratio i/j = C2/C1. By defining the auxiliary variable dy = cidt, the rate equation for the 
concentration ratio is 

^ = -{D2-D,){^ + r). (3) 

Because i/j is a strictly decreasing function of y for D2 > -Di , the asymptotic behavior may 
be estimated by ignoring the quadratic term in Eq. (3). The corresponding solution is 



simply 



~ e 



-{D2-D1) J*dt' Cl(t') 



(4) 



Two possibilities for the behavior of y = J dt' Ci{t') as t 00 can occur. If y reaches 
a finite limit as t ^ 00, then C2{t) cx ci{t); further, ci{t) and C2{t) must decay as t~^. 
It is immediate to show that these two conditions are incompatible with the initial rate 
equations. On the other hand, if y — > 00 as t — > 00, then C2/C1 — > 0, asymptotically. Thus 
in Eq. (2a), C1C2 may be neglected compared to cf and the resulting approximation gives 
ci{t) ~ {2Dit)~^. Using this in the equation for ij;, one ultimately finds 



C2 ~ C(2Dit)-(i+^), C 



C2(0) 

ci(0) 



1 + C2(0)/Ci(0) 

ci(0) 



-I 5 



(5) 



for t — > 00, with 5 = {D2 — Di)/2Di. While the less mobile species decays as in the 
long-time limit, as in homogeneous single-species annihilation, the more mobile species 
decays non-universally as with a = 1 + 6 = {Di + D2)/2Di. 

It is worth noting that when S is small, the asymptotic behavior given in Eq. (5) 
is reached only at very long times. To demonstrate this, consider, e. g., equal initial 
concentrations of the two species. Then solving Eq. (3) with ip{0) ~ 1 and < 5 <C 1 gives 
i^iv) — C2/C1 ~ 1 + 51n[ci(|/)/ci(0)]. Substituting this into Eqs. (2), one finds 



Cl,2{t) 



4Dit 



l±-ln{D,t) 



(6) 



for 1 <^ i <^ e^/*^. Consequently, the crossover time between the intermediate asymptotics, 
Eq. (6), and the final asymptotics, Eq. (5), diverges as e^/"^ for 6^0. 

For a finite number of species, the kinetics is similar to that of the two-species case. 
Namely the concentration of the least mobile species decays as while the more mobile 
species exhibit non- universal power-law decays with associated exponents {Di + Dk)/2Di, 
where Dk is the diffusion coefficient of the k^^ species. 

When the diffusivities of the species are drawn from a continuous distribution, the 
rate equation for P{D,t), the concentration of species with diffusivity D, becomes the 
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integro-difFerential equation 



oo 



dP{D,t) 
di 



P{D, t) I dD' {D + D') P{D', t) 



(7) 







P{D,t)[DVo{t)+Vi{t)]. 



Here Vk{t) denotes the k^^ moment of the difFusivity distribution P{D,t), 



oo 




dDD''P{D,t). 



(8) 







Note that the zeroth moment of the diffusivity distribution is just the particle density, 
c{t) = Vo{t): while the average diffusion coefficient is expressed in terms of the zeroth and 
first moment by {D) = Vi{t)/Vo{t). 

Equations similar to Eq. (7) describe the kinetics of several irreversible processes, such 
as diffusion- reaction aggregation [7] , ballistic annihilation [4] , and ballistic aggregation [8] . 
For the aggregation problem, a reaction rate of the form {D + D') is known as the "sum" 
kernel, for which the rate equations are exactly soluble [7] . In heterogeneous annihilation, 
the rate equations are also soluble by elementary analysis. A formal, but implicit solution 
to Eq. (7) is 



"-0 

Further, by integrating Eq. (7) over D, the following differential equation which relates 
the moments ■Po(^) and Vi{t) is obtained 



With the initial condition Po(^ = 0) = 1; which fixes the scale of the initial distribution, 
the solution to Eq. (10) is 





dVo 



= -2VoVi. 



(10) 



dt 








(11) 
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Combining Eqs. (9) and (11) leads to a simplification in which only the unknown zeroth 
moment, or particle density c{t) — 'Po{t), appears in the formal solution 

t 



P{D, t) = P{D, 0) Vc(t) exp 



-D J dt' c{t') 







(12) 



A scaling analysis of this solution indicates that two types of behavior can occur: 
one for initial diffusivity distributions with a finite non-zero lower cutoff, Dmin > 0, and 
the other for situations where -Dmin = 0. For both cases, we consider only those initial 
distributions which are homogeneous near -Dmin, namely, P{D, 0) ~ (D — -Dmin)'^ as -D — > 
-Dmin with > — 1 for normalizability. This restriction leads to mathematical tractability 
as well as being a natural illustrative choice. 

For Dmin = 0, the average diffusion coefl&cient clearly decay to zero, and under mild 
restrictions, the particle concentration also decays to zero. For an initial distribution of 
the form P{D, 0) ~ D^^ as D ^ 0, it is natural to assume power-law decays in the average 
concentration and average diffusivity: c ~ and (D) ~ for t — > oo. Correspondingly, 
the time-dependent diffusivity distribution is expected to approach the scaling form, 

P{D,t)c^t^-°'^Dt^). (13) 

Upon substituting this scaling ansatz into Eq. (12), one finds two relations which then 
determine the exponents a and /3. First, for the argument of the exponential in Eq. (12) 
to be dimensionless, if follows that a + (3 — 1. Similarly, for the prefactor to be a function 
only of the scaling variable x = Dt^, the relation a = 2(3{1 + /j,) must hold. These two 
conditions determine the dependence of the fundamental exponents on fi , 

2 + ^ 1 ,,,, 

a= 3^ , (14) 

while the scaling function is 

^x) = x^e-"". (15) 

Let us now justify the existence of the scaling ansatz by constructing explicit solutions 
to Eq. (12). Consider the initial distribution P(-D, 0) = D^'e-^/V{^i+ 1), where the 



numerical factor normalizes the initial density to unity, for convenience. Substituting this 
initial distribution into Eq. (12) and then integrating the resulting equation over D, one 
finds the simple differential equation 

du 



dt 



u 



(16) 



for the auxiliary function u{t) — 1 -\- j dt'c{t'). With the initial condition, u(t = 0) = 1, 



Eq. (16) is readily solved, from which the concentration is 



, . du \ 
c{t) = -^= [l + (2// + 3)t 



3+2^ 



(17) 



Substituting this in Eq. (12), we obtain for the diffusivity distribution 



P{D,t) 



1 + (2// + 3)t 



l + M 

3 + 2m 



exp —D 



1 + (2// + 3)t 



3+2/j 



(18) 



r(/x+i) 

As expected, this solution approaches the scaling form of Eq. (13) asymptotically, with the 
exponents and scaling function given by Eqs. (14) and (15), respectively. 

Consider now initial distributions with Dmin > 0. For concreteness, we examine the 
case where 

After substitution of this initial distribution in Eq. (12), the analog of Eq. (16) is 

du 



dt 



^-2M-2g-2r>^in('U-l) _ 



(20) 



Solving this equation in the limit of t — > oo, the asymptotic form for the density is 

1 



c{t) 



2D • / 



^ 2 + 2/x 



(21) 



while the full distribution of diffusivities has the form 



P{D,t) 



2D ■ 



(22) 



p(AO) 

Thus the concentration of the species with diffusivity D decays algebraically with a charac- 
teristic non-universal exponent {D + Dmin)/2Dmin5 but also with a logarithmic prefactor 



This distribution can be expressed in a scaling form by identifying the appropriate 
scahng variable. For this identification, it is helpful to consider the average diffusion 
coefficient which has the time dependence (D) — Dmin ~ -Dmin/ ln(2Dmint). This suggests 
that the scaling variable is 

.= i^#^P-D„..). (23) 

With this identification and in the scaling limit t ^ oo and D £>min but with x finite, 
the asymptotic distribution Eq. (22) can be written in the form 

(24) 

with scaling function ^{x) again equal to x^^e~^. Thus for both -Dmin = and -Dmin > 0, 
the diffusivity distribution approaches the scaling form, given by Eq. (13) and Eq. (24), 
respectively, with characteristic width ^"1/(2+^/^) in the former case and l/ln(t) in the 
latter. 

III. Heterogeneous Annihilation in Low Dimensions 

In the diffusion-controlled limit, fluctuations effects govern the long-time kinetics of 
heterogeneous annihilation for (i < 2, in close analogy with homogeneous single-species 
annihilation [1]. A method which is ideally suited to account for the kinetics in this 
regime is the Smoluchowski theory [5]. We adapt this approach to treat the kinetics of 
heterogeneous annihilation when d < 2. To describe the method, consider first the simple 
and illustrative example of a two-species system in one dimension in which a background 
of identical particles with diffusivity D and density c contains relatively rare impurities of 
diffusivity Dj at concentration c/ <^ c. 

In the Smoluchowski approach, we first compute particle flux to a "reference" ab- 
sorbing particle due to the rest of the particles which comprise the uniform background. 



The requisite solution to the background concentration is c(x,t) = Cooerf (a:/v4Dt), from 
which the particle flux at the reference particle is = Cqo \/d7^- This is identified as the 
effective microscopic reaction rate, k. 
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For the two-component system of background and impurity, there are distinct rates 
associated with reactions between background particles and between the background and 
impurities. In the hmit of low impurity density, we neglect the influence of background- 
impurity reactions on the background density, as well as reactions among the impurities. 
Under these restrictions, the rate equations become 



where kBB and ksi are the effective rates for background-background and background- 
impurity reactions. Note that the numerator inside the square root involves the relative 
diffusivities of the two reacting species. From the first equation, the background concen- 
tration vanishes as 



This deviates from the exact result [2, 3] for single-species annihilation by the numerical 
factor c(t)/cexact(^) = 7r/2. The crucial feature of the Smoluchowski approach is that by 
the form of the rate equation for c/, the coefficient of c{t) determines the exponent of the 
decay of the impurity species. We thereby find 



with e = Di/D. 

As we discuss in the next section, this prediction agrees rather well with numerical 
results. Intuition for this non-universal behavior can be gained by considering the impurity 
problem in the case where particles coalesce, rather than annihilate. This situation is con- 
siderably simpler than that of annihilation, as the enclosing "cage" , defined by the nearest 
neighbors of the impurity, evolves only by diffusion. Thus the many-body problem may 
be reduced to the three-body problem of the impurity and its two nearest-neighbors. By 
exploiting simple geometric equivalences, this problem can be transformed to the survival 
of a single random walker which diffuses within an absorbing two-dimensional wedge whose 




(25) 



c{t) = y/7r/32Dt. 



(26a) 



c/(t)~t-V(i+^)/8, 



(266) 
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opening angle depends on D/Dj. In this latter problem, the survival probability decays 
as t~'^, with a = 7r/(2 cos~^(e/(l + e))). Here 9 = cos"-*^ [e/(l + e)] is the opening angle of 
the wedge [9]. 

Such a rigorous mapping does not exist for the impurity problem when particles 
annihilate. This is an intrinsically a many-body process, since the cage surrounding a 
given particle can involve distant neighbors. Nevertheless, the Smoluchowski approach 
is essentially identical for both the annihilation and the coalescence reactions, except for 
overall factors of 2 in the rate equations. Thus the equivalent of Eq. (26b) for aggregation 
is ciit) ~ t While this exponent value deviates considerably from the exact 

result give above, the e dependence of the exponent is qualitatively correct. 

The mechanism underlying the non-universal decay of the impurity is the equivalence 
to the survival probability of a diffusing particle inside an absorbing interval whose length 
L grows as . Since the probability density of the particle spreads over a spatial extent 
of the order of \fDi, the survival probability decays exponentially for /3 < 1/2. However, 



in the marginal case (3 = 1/2, i. e., L ~ Ay/Dt, the survival probability S{t) decays 




mensionless parameter A. In the reaction process, the size of the cage which surrounds 



are important because they determine the amplitude A in the growth of the cage. As 
shown above, the phenomenological Smoluchowski treats annihilation and coalescence on 
the same footing and thus provides a convenient description of the non-universal behavior. 

We have also explored the case of small difference between diffusion coefficients, < 
5 <^ 1. We find that intermediate asymptotic behavior initially occurs before the final 
asymptotic of Eq. (26) sets in. In particular, for initially equal concentrations 




algebraically in time, S{t) ~ t 



but with the decay exponent dependent on the di 



the impurity also grows as t^/^ for both annihilation and coalescence, since the overall 
density is decaying as t~^/'^. However, the microscopic differences in the two reactions 




(27) 



for 1 <^ t -C e^/*^. The crossover time between intermediate asymptotics Eq. (27) and final 
asymptotics Eq. (26) diverges as e^/*^ for (5 — > 0. 
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Let us now consider heterogeneous annihilation in one dimension with a continuous 
distribution of particle difFusivities. Since the Smoluchowski approach in one-dimension 



yields a flux that varies as oc JDjt, one finds a rate equation of the form of Eq. (7) 



with the kernel K(D,D') oc ■\/{D'^\^W)Jt- The correspondence can be made even closer 



by introducing the modified time variable T = Ay^t/TT which eliminates the explicit time 
dependence. The rate equation for P{D,T) becomes 



dP{D,T) 

df 



-P{D,T) J dDWD + D' P{D',T). (28) 



We are unable to find either a scaling solution to this equation or the exponents. 
However, presumably exact values of the exponents can be obtained by consideration of a 



closely related and more tractable model. We replace the kernel K{D, D') = \/D + D' by 
one with the same homogeneity degree, K{D, D') = \fD + \fD' . The two kernels obey the 
bounds ^yl^ < \/D -\- D' < ^/D + a/D'. This suggests that the asymptotic behavior 
from the initial kernel should be identical to that predicted by the simpler sum-root kernel 
K{D,D') = \fD^^. 

With this modified kernel, the rate equation becomes 

^^^f^ = -P{p,T) [x/^Po(T) +Pi/2(T)], (29) 

oo 

with Vi/2{T) = / dD y/D P{D,T). Repeating the steps employed for the mean-field 



treatment, the formal exact solution to Eq. (29) is 



P(A r) = P(A 0) V^exp 



1 

-D J dT'c{T') 



(30) 



Again we substitute the scaling ansatz P{D, t) ~ t^~°'^{Dt^) into the above rate equation 
with the initial condition P(D, 0) ~ D^^ as D ^ 0. By expressing all factors in a scaling 
form we find the two exponent relations 2a + P = 1 and a = 2^(1 + n) and hence, 

« = 7 — ^ = 7 — • (31) 
5 + 4/i 5 + 4/i ^ ^ 
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Here the exponents a. and /? refer to the time dependence of c(t) and (D) in terms of the 
physical time variable t. 

These results can be straightforwardly generalized to arbitrary d < 2 within the 
Smoluchowski approach. A simple calculation shows that the flux to a trap behaves as 
£)d/2^-i+d/2 £qj. < 2. In the rate equations, this suggests the introduction of the modified 
time T ~ (Dt)^/'^. One thereby obtains an equation similar to Eq. (28), but with the 
kernel K{D,D') = (D + D')'^^'^. Next we replace this kernel by more tractable form, 
K[D,D') = D*^/^ + D"^/^. Applying a scaling analysis to the rate equations with this 
kernel then gives the exponents a = 2d{l + iJ,)/[d + 4(1 + //)] and P = d/[d + 4(1 + //)]. 

In two dimensions, the Smoluchowski approach shows that the reaction kernel is 
K{D, D') — At:{D + D') / \n.{Dt). In this sense, two dimensions is the marginal case which 
demarcates the regime where the reaction rate is time independent (for d > 2) from the 
regime where the reaction rate varies as a power law in time {d <2). The corresponding 
rate equation for P(-D, t) in d = 2 is 



We analyze this equation in the same spirit as that employed for d 7^ 2. First, we approx- 
imate the slowly varying logarithmic factor by ln[(D + D')t\ ~ ln[(Z))t] ~ aln(t) and then 
introduce the auxiliary time variable T ATrt/alii{t), so that Eq. (32) becomes identical 
in form to the mean- field rate equation, Eq. (7). Thus replacing t by T, the results of 
the preceding section apply to the present situation. In particular, at asymptotically large 
times P{D, t) approaches the scaling form 



with exponents and the scaling function given by Eqs. (14) and (15). 

IV. Numerical Results 

Consider the specific and illustrative example of a two-component system consisting 
of an impurity particle with diffusivity Dj within a background of identical particles with 



00 




(32) 







(33) 
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difFusivity D. We first describe the results of a time series analysis for the survival prob- 
ability of a static impurity. Then we discuss complementary Monte Carlo simulations for 
both the time dependence of the survival probability and the eventual survival probability 
in a finite system. Both techniques yield estimates for the decay exponent of the impurity 
that are in good agreement with the Smoluchowski theory predictions. 

The time series expansion is a general technique for evaluating numerically the initial 
terms in the exact power-series expansion in time for various observables that characterize 
reactive systems [10]. The technique involves the successive application of a suitably- 
defined evolution operator on the initial state of the system. After each application of 
the evolution operator, configurations for which the impurity particle has not reacted 
contribute to the survival probability. To compute the first n coefficients of the time power 
series, it is sufficient to consider the evolution of an n-site ring. The primary limitation of 
the technique is computer memory, since the number of configurations grows exponentially 
with the order of the expansion. For the case of the static impurity, we obtained the 
expansion to order 21, as shown in Table 1. We then apply the Pade analysis method 
suggested in Ref. [11] to estimate the asymptotic properties of the survival probability. In 
Figure 1, the diagonal Pade approximants for the survival probability, S[n,n]it) are plotted 
forn = 8,9, and 10. These approximants are essentially identical for t < 6 and we conclude 
that they accurately describe the true survival probability within this time range. These 
Pade approximants are then fitted to the power-law form S{t) ~ t~" using the fitting 
procedure also suggested in [11]. We thereby estimate the exponent value a = 0.38 ± 0.01 
(Table 2). 

Monte Carlo simulations were also performed for the impurity problem in one dimen- 
sion. A typical measurement involved 10 configurations of a periodic chain of 5 x 10^ sites 
with initial concentrations c/(0) = 0.01 and c(0) = 0.99. The data for the time depen- 
dence of both ci{t) and c{t) arc quite straight on a double logarithmic plot and exponent 
estimates may be made based on the local slopes of nearby data points. From the slopes of 
first-neighbor, second-neighbor, and third-neighbor pairs, we estimate that the background 
species decays as t~°'^ , where as = 0.50 ± 0.01. Here the error bar indicates the magni- 
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tude of the fluctuations in the local slopes in the time regime where the data is most linear 
(between 20 ^ t <^ 10,000). Since ckb = 1/2, this analysis provides a useful confidence 
test. For the impurity species, the same analysis method suggests a — 0.37 ±0.01 (Fig. 2). 

A complementary approach is based on a finite size scaling analysis for the eventual 
survival probability of an impurity, S{L) = \imt^^ S{L,t), when starting with a single 
impurity and the rest of the sites occupied by background particles on finite ring with an 
odd number of sites, L. We postulate that S{L) ~ (clearly 7 = 1 for Dj — D), while 
for short times, the survival probability must be independent of the size of the system, so 
that S{L,t) ~ t~°'. Since the two basic length scales are y/Dt and L, we assume that the 



survival probability has the scaling form S{L,t) ~ f{-\/Di/ L). The scaling function / 

must have the limiting behaviors /(x) ~ 1 for x ^ and /(x) ~ x^*^ for x ^ 00 to account 
for the two asymptotic limits. These conditions imply that 7 = 2a. An advantage of the 
measurement of the eventual survival probability is that it typically gives more accurate 
exponent estimates compared to the time-dependent observations. 

We have performed Monte Carlo simulations of finite single impurity systems until 
the impurity disappears or a single impurity remains to measure the exponent 7 directly. 
In Figure 3, we present the average survival probability for 10^ realizations, for systems of 
size 2*^ + 1, with n = 1, . . . , 8. The data suggests the decay exponent of 7 = 0.75 ±0.01, 
consistent with the series and Monte Carlo estimates of a = 0.38 ±0.01. Table 2 compares 
the decay exponents found for the other values of e = Dj /D to the corresponding exponent 
given by the Smoluchowski approach. The Smoluchowski approximation yield a remarkably 
good description of the asymptotic behavior of the impurity decay in the annihilation 
reaction. 



We have shown that there is non-universal kinetic behavior in diffusion-controlled 
heterogeneous single-species annihilation. Typically, the concentration of the species with 
the smallest diffusion coefficient has the same time dependence as in homogeneous single- 
species annihilation, ^ ± ^ ^ 0. However, the concentrations of the more mobile species 




V. Discussion 
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decay as power laws in time, with non-universal exponents that depend on the ratios of 
the corresponding diffusivities to that of the least mobile species. These exponents were 
determined by a mean-field approximation, which should be valid for d > 2, and by a 
phenomenological Smoluchowski theory which should apply for (i < 2. 

Our numerical studies focused on the limiting case where there is an infinitesimal 
concentration of one species in a homogeneous background of another species. For this 
impurity problem, the behavior of the cage which surrounds the impurity is the crucial 
ingredient which determines the kinetics of the impurity. For the coalescence reaction, 
A + A^A, the cage consists of the two nearest-neighbors of the impurity. These enclosing 
particles continue to undergo diffusion even if there are reactions with other more distant 
particles in the system. This allows one to reduce the many-body system to a three-particle 
system which can then be solved exactly in one dimension. In the annihilation reaction, a 
nearest-neighbor of the impurity can disappear because of a reaction with the next-nearest 
background particle. This can cause a relatively large rearrangement in which the cage 
expands to the next background particle. 

The Smoluchowski approach turns out to be ideally suited for providing a simple 
description of the non-universal behavior of the impurity survival probability. Although 
this approach involves an uncontrolled approximation, the Smoluchowski method provides 
a simple way to quantify the enclosing cage of the impurity for both the coalescence and 
annihilation reactions. This is a remarkable success in view of the relatively complex cage 
dynamics in the annihilation reaction. From the information about the size of the cage, the 
decay exponent follows directly. In particular, for the annihilation reaction, the exponent 
predictions compare well with numerical estimates. 
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Figure Captions 

Figure 1. Time dependence of the diagonal Fade approximants S[n,n]{t) = 8,9,10) 
based on the time series for the survival probability in the case of a static 
impurity. A line of slope —0.38 is plotted for reference. 

Figure 2. Monte Carlo simulation results for the survival probability of impurity (squares) 
and background particles (circles) based on 10 configurations of a one-dimensional 
chain of 500,000 sites. The initial impurity and background concentrations are 
0.01 and 0.99, respectively. 

Figure 3. The eventual survival probability S{L) (circles) of a single static impurity on 
a one-dimensional ring of L sites in which the rest of the ring is initially filled 
with diffusing particles. The data shown is based on 10^ realizations. A line of 
slope —0.75 is shown for reference. 

Table Captions 

Table 1. The first 21 coefficients in the time series for the survival probability of a static 

impurity, defined by S{t) = c^t". 
Table 2. Comparison of the numerical estimates for the exponent 7 that characterizes 

the L dependence of the eventual survival probability with the corresponding 

predictions from the Smoluchowski theory. Estimates are based on averaging 

over 10^ realizations. 



18 



